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We investigate diffusion in supersonic, turbulent, compressible flows. Supersonic turbulence can be 
characterized as network of interacting shocks. We consider flows with different rms Mach numbers 
and where energy necessary to maintain dynamical equilibrium is inserted at different spatial scales. 
We find that turbulent transport exhibits super-diffusive behavior due to induced bulk motions. 
In a comoving reference frame, however, diffusion behaves normal and can be described by mixing 
length theory extended into the supersonic regime. (Accepted for publication in Phys. Rev. E) 
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*0 \ I. INTRODUCTION 

Laboratory and terrestrial gases and liquids are usually well described by incompressible flows (e.g. |2fij|L In contrast, 
' the dynamical behavior of typical astrophysical gases, are characterized by poorly understood highly compressible 
supersonic turbulent motion (see e.g. |58j). For example, the large observed lincwidths in large molecular clouds show 
direct evidence for the presence of chaotically oriented velocity fields with magnitudes in excess of the sound speed. 
This random motion carries enough kinetic energy to counterbalance and sometimes overcompensate the effects of 
self-gravity of these clouds |60j. The intricate interplay between supersonic turbulence and self-gravity determines 
\ the overall dynamical evolution of these clouds and their observable features such as their density structure, the star 
. formation rate within them, and their lifetimes. Thus, it is importance for the description of many astrophysical 
systems to understand in detail the momentum and heat transfer properties of compressible turbulent gases. 

Some important clues on the nature and efficiency of mixing associated with the clouds' supersonic turbulence can 
be constrained by the observed metallicity distribution of the stars formed within them. In the Pleiades cluster, stars 
which emerged from the same molecular cloud have nearly identical metal abundance |59j | . This astronomical context 
(-^ ■ therefore imposes a strong motivation for a general analysis of the transport and mixing processes in compressible 
^ | supersonically turbulent media. 

Analytical and numerical studies of diffusion processes are typically restricted to certain families of statistical 
processes, like random walk |33| or remapping models or certain Hamiltonian syste ms Il9l. The direct numerical 
modeling of turbulent physical flows mostly concentrates on incompressible media (e.g. |!4ll35l l38| ). but some studies 
have been extended into the weakly compressible regime [Tj, [l^, IS H2J Eg • Although highly c omp ressible supersonic 
turbulent flows have been studied in several specific astrophysical contexts 0, H |1 H [ill 111 13 13 HI El IH E 
Ell E3 . liM Eil I45L I52I I53I IHil l56l IH3 , the diffusion properties of such flows have not been investigated in detail. 

It is the goal of this paper to analyze transport phenomena in supersonic compressible turbulent flows and to 
demonstrate that - analogous to the incompressible case - a simple mixing length description can be found even for 
strongly supersonic and highly compressible turbulence. We first briefly recapitulate in Section[n]the Taylor formalism 
for describing the efficiency of turbulent diffusion in subsonic flows. In Section IlIII we describe the numerical method 
which we use to integrate the Navier-Stokes equation. In Section Hvl we report the diffusion coefficient obtained in our 
numerical models, and in Section[V]we introduce an extension of the well known mixing length approach to diffusion 
into the supersonic compressible regime. Finally, in Section IVII we summarize our results. 
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II. A STATISTICAL DESCRIPTION OF TURBULENT DIFFUSION 



Transport properties in fluids and gases can be characterized by studying the time evolution of the second central 
moment of some representative fluid-elements' displacement in the medium, 

g(t-t') = {[r i (t)-?i(t')} 2 )i, (1) 

where the average (•).; is taken over an ensemble of passively advected tracer particles i (e.g. dye in a fluid, or smoke 
in air) that are placed in the medium at a time t' at positions ri(t'); or where the average is taken over the fluid 
molecules themselves (or equivalently, over sufficiently small and distinguishable fluid elements). The dispersion in 
one spatial direction, say along the ^-coordinate, is £%(t — t') = {[xi(t) — Xi(t')] 2 )i. For isotropic turbulence it follows 
that Q. — & = = 1/3 £2. For fully-developed stationary turbulence, the initial time t' can be chosen at random 
and for simplicity is set to zero in what follows. 

The quantity £p(t) can be associated with the diffusion coefficient D as derived for the classical diffusion equation, 

where n(ri,t) is the probability distribution function (pdf) for finding a particle i at position fi(t) at time t when it 
initially was at a location r*i(0). This holds if the particle position is a random variable with a Gaussian distribution 
In the classical sense, n{f,t) may correspond to the contaminant density in the medium. Equation [2] holds for 
normal diffusion processes and for time scales larger than the typical particles' correlation time scale t. 
In general, however, the Lagrangian diffusion coefficient is time dependent and can be defined as 

D(t) = &-=2(f i (t).v i (t))i, (3) 

where = dfi(t)/dt is the Lagrangian velocity of the particle. The diffusion coefficient along one spatial direction, 
say along the x-coordinate, follows accordingly as D x = d^(t)/dt = 2(xi(t)v x .i(t))i. Equation^Jholds for homogeneous 
turbulence with zero mean velocity. From ri(i) = f*j(0) + L Vi{t')dt' it follows that 



D(t) = 2 



;(0)+ / Vi(t')dt' 



Vi(t)) = 2 (viit 1 ) ■ Vi(t))idt' . (4) 







The above expression allows us to related D(t) to the trace of the Lagrangian velocity autocorrelation tensor trC(£ — 
t') = {vi(t') .€4(t))i as 

D(t) = 2 1 tvC{t- t')dt' =2 1 txC{t')dt' , (5) 
Jo Jo 

a result which was derived by Taylor already 1921 [5{|. This formulation has the advantage that it is fully general 
and that it allows us to study anomalous diffusion processes. Note, that strictly speaking any transport process with 
&{t) not growing linearly in time is called anomalous diffusion. This is always the case for time intervals shorter 
than the correlation time r, but sometimes anomalous diffusion can also occur for t 3> r. If £p(t) oc t a and if a < 1 
transport processes are called subdiffusive, if a > 1 they are called superdiffusive [T2I ITiil l2ri |27| . Studying transport 
processes directly in terms of the particle displacement, i.e. Equation ^ is useful when attempting to find simple 
approximations to the diffusion coefficient D(t) for example in a mixing length approach. 



III. NUMERICAL METHOD 



In order to utilize the above formalism, we carry out a series of numerical simulation of supersonic turbulent flows. 
A variety of numerical schemes can be used to describe the time evolution of gases and fluids. By far the most 
widely-used and thorou ghly -studied class of methods is based on the finite difference representations of the equations 
of hydrodynamics (e.g. 49j ) . In the most simple implementation, the fluid properties are calculated on equidistant 
spatially fixed grid points in a Cartesian coordinate system. Finite difference schemes have well defined mathematical 
convergence properties, and can be generalized to very complex, time varying, non-equidistant meshes with arbitrary 
geometrical properties. However, it is very difficult to obtain a Lagrangian description, which is essential when dealing 
with compressible supersonic turbulence with a high degree of vorticity. Methods that do not rely on any kind of 
mesh representation at all are therefore highly desirable. 
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For the current investigation we use smoothed particle hydrodynamics (SPH), which is a fully Lagrangian, particle- 
based method to solve the equations of hydrodynamics. The fluid is represented by an ensemble of particles, where flow 
properties and thermodynamic observables are obtained as local averages from a kernel smoothing procedure (typically 
based on cubic spline functions) fft. l34|. Each particle i is characterized by mass TOj, velocity Vi and position and 
carries in addition density pi 1 internal energy ej or temperature Tj, and pressure pi. The SPH method is commonly 
used in the astrophysics community because it can resolve large density contrasts simply by increasing the particle 
concentration in regions where it is needed. This versatility is important for handling compressible turbulent flows 
where density fluctuations will occur at random places and random times. The same scheme that allows for high spatial 
resolution in high-density regions, however, delivers only limited spatial resolution in low-density regions. There, the 
number density of SPH particles is small and thus the volume necessary to obtain a meaningful local average tends 
to be large. Furthermore, SPH requires the introduction of a von Neumann Richtmyer artificial viscosity to prevent 
interparticle penetration, shock fronts are thus smeared out over two to three local smoothing lengths. Altogether, 
the performance and convergence properties of the method are well understood and tested against analytic models 
and other numerical schemes, for example in the context of turbulent supersonic astrophysical flows |2lL 123. 12^. l8ll| . 
and its intrinsic diffusivity is sufficiently low to allow for the current investigation of turbulent diffusion phenomena 

M 

To simplify the analysis we assume the medium is infinite and isotropic on large scales, and consider a cubic volume 
which is subject to periodic boundary conditions. The medium is described as an ideal gas with an isothermal equation 
of state, i.e. pressure p relates to the density p as p — c^p with c s being the speed of sound. Throughout this paper we 
adopt normalized units, where all physical constants (like the gas constant), total mass M, mean density (p), and the 
linear size L of the cube all are set to unity. The speed of sound is c s = 0.05, hence, the sound crossing time through 
the cube follows as £ SOU nd = 20. In all models discussed here, the fluid is represented by an ensemble of 205 379 SPH 
particles which gives sufficient resolution for the purpose of the current analysis. 

Supersonic turbulence is known to decay rapidly 0, IjnJ |3(| El Ell • Stationary turbulence in the interstellar 
medium therefore requires a continuous energy input. To generate and maintain the turbulent flow we introduce 
random Gaussian forcing fields in a narrow range of wavenumbers such that the total kinetic energy contained in 
the system remains approximately constant. We generate the forcing field for each direction separately and simply 
add up the three contributions. Thus, we excite both, solenoidal as well as compressible modes at the same time. 
The typical ratio between the solenoidal and compressible energy component is between 2:1 and 3:1 in the resulting 
turbulent flow (see e.g. Figure 8 in |23|.). We keep the forcing field fixed in space, but adjust its amplitude in order 
to maintain a constant energy input rate into the system compensating for the energy loss due to dissipation (for 
further details on the method see |H,|2^|). This non-local driving scheme allows us to exactly control the (spatial) 
scale which carries the peak of the turbulent kinetic energy. It is this property that motivated our choice of random 
Gaussien fields as driving source. In reality the forcing of turbulence in the interstellar medium is likely to be a 
multi-scale phenomenon with appreciable contributions from differential rotation (i.e. shear) in the Galctic disk and 
energy input from supernovae explosions ending the lives of massive stars 1^3 . Comparable to the values observed 
in interstellar gas, we study flows with Mach numbers in the range 0.5 to 10, where we define the Mach number 
from the one- dimensional rms velocity dispersion o~ v as M. = o~ v /c s . For each value of the Mach number we consider 
three different cases, one case where turbulence is driven on large scales only (i.e. with wavenumbers k in the interval 
1 < k < 2), intermediate- wavelength turbulence (3 < k < 4), and small-scale turbulence (7 < k < 8), as summarized 
in TableQ] Note that our models are not subject to global shear because of the adopted periodic boundary conditions. 
We call turbulence "large scale" when the Fourier decomposition of the velocity field is dominated by the largest 
scales possible for the considered volume L 3 , i.e. the system becomes isotropic and homogeneous only on scales larger 
than L. On scales below L it may exhibit a considerable degree of anisotropy. This is most noticeable in the case 
I < k < 2, because wavenumber space is very poorly sampled and variance effects become significant. The system is 
dominated by one or two large shock fronts that cross through the medium. In the interval 7 < k < 8 the number 
of Fourier modes contribution to the velocity field is large, and the system appears more isotropic and homogeneous 
already on distances smaller than L. This trend is clearly visible in Figure 

Similar to any other numerical calculations, the models discussed here fall short of describing real gases in com- 
prehensive details as they cannot include all physical processes that may act on the medium. In interstellar gas 
clouds, transport properties and chemical mixing will not only be determined by the compressible turbulence alone, 
but the density and velocity structure is also influenced by magnetic fields, chemical reactions, and radiation transfer 
processes. Furthermore, all numerical models are resolution limited. The turbulent inertial range in our large-scale 
turbulence simulations spans over about 1.5 decades in wavenumber. This range is considerably less than what is 
observed in interstellar gas clouds. The same limitation holds for the Reynolds numbers achieved in the models, 
they fall short of the values in real gas clouds by several orders of magnitude. Nevertheless, despite these obvious 
shortcomings, the results derived here do characterize global transport properties in interstellar gas clouds and in 
other supersonically turbulent compressible flows. 
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TABLE I: Model properties 
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(2) 


(3) 


(4) 
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(9) 


(10) 


(11) 


(12) 


(13) 


model 
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M 






(J y 






Dy (OO) 


D z {oo) 




/k 


2S y /k 


2a z /k 


01 


1. 


.2 


0.6 


35.3 


0.030 


0.028 


0.027 


0.027 


0.021 


0.019 


0.030 - 


0.060 


0.028 - 0.057 


0.027 - 


0.054 


Oi 


3. 


.4 


0.5 


39.1 


0.026 


0.026 


0.025 


0.010 


0.010 


0.009 


0.013 - 


0.017 


0.013 - 0.017 


0.013 - 


0.017 


Os 


7. 


.8 


0.4 


46.2 


0.021 


0.022 


0.022 


0.005 


0.005 


0.005 


0.005 - 


0.006 


0.005 - 0.006 


0.005 - 


0.006 


1£ 


1. 


.2 


1.9 


10.4 


0. 106 


0.084 


0.098 


0. 140 


0.069 


0.111 


0.106 - 


0.213 


084 — 1 67 


0.098 - 


0.196 


li 


3. 


.4 


1.9 


10.6 


0.097 


0.096 


0.092 


0.042 


0.047 


0.038 


0.048 - 


0.065 


0.048 - 0.064 


0.046 - 


0.061 


Is 


7. 


.8 


1.7 


11.5 


0.086 


0.089 


0.087 


0.025 


0.026 


0.024 


0.021 - 


0.024 


0.022 - 0.025 


0.022 - 


0.025 


21 


1. 


.2 


3.1 


6.5 


0.173 


0.129 


0.158 


0.223 


0.103 


0.169 


0.173 - 


0.346 


0.129 - 0.257 


0.158 - 


0.315 


Si 


3. 


.4 


3.1 


6.4 


0.167 


0.155 


0.151 


0.084 


0.071 


0.063 


0.083 - 


0.111 


0.077 - 0.103 


0.075 - 


0.100 


2s 


7. 


.8 


3.2 


6.3 


0.154 


0.163 


0.157 


0.044 


0.054 


0.047 


0.038 - 


0.044 


0.041 - 0.046 


0.039 - 


0.045 


31 


1. 


.2 


5.2 


3.8 


0.301 


0.252 


0.227 


0.314 


0.245 


0.169 


0.301 - 


0.603 


0.252 - 0.505 


0.227 - 


0.454 


Si 


3. 


.4 


5.8 


3.5 


0.261 


0.287 


0.316 


0.131 


0.189 


0.233 


0.130 - 


0.174 


0.143 - 0.191 


0.158 - 


0.211 


3s 


7. 


.8 


5.8 


3.4 


0.297 


0.288 


0.289 


0.106 


0.092 


0.091 


0.074 - 


0.085 


0.072 - 0.082 


0.072 - 


0.083 



4£ 1..2 8.2 2.4 0.467 0.318 0.444 0.693 0.241 0.558 0.467 - 0.933 0.318 - 0.635 0.444 - 0.887 
4% 3. .4 9.7 2.1 0.451 0.478 0.520 0.248 0.323 0.349 0.225 - 0.301 0.239 - 0.319 0.260 - 0.347 
4s 7. .8 10.4 1.9 0.532 0.513 0.519 0.194 0.167 0.170 0.133 - 0.152 0.128 - 0.147 0.130 - 0.148 



1. column: Model identifier, with the letters £, i, and s standing for large-scale, intermediate-wavelength, and short-wavelength turbulence, 

respectively. 
2. column: Driving wavelength intervall. 
3. column: Mean Mach number, defined as ratio between the time-averaged one-dimensional velocity dispersion a v — 3~ 1 ^ 2 (a"^ + tx^ + ct^) 1 ' 2 
and the isothermal sound speed c s , A4 — a v /c s . The valu es fo r the different velocity components x, y, and z may differ considerably, especially 
for large-wavelength turbulence. Please recall from Section lllll that the speed of sound is c s — 0.05, and thus the sound crossing time i souri d — 20. 

4. column: Average shock crossing time through the computational volume. 
5. to 7. column: Time averaged velocity dispersion along the three principal axes x, y, and z, e.g. for the x-component 

9l = f*{(v*Ai')-(v*.i(t'))i) a )idt'/t. 
8. to 10. column: Mean-motion corrected diffusion coefficients along the three principal axes computed from Equation[3]ior time intervals t ;§> t. 
11. to 13. column: Predicted values of the mean motion corrected diffusion coefficients D x , D , and D z from extending mixing length theory into 

the supersonic regime (SectionQ. 

IV. TRANSPORT PROPERTIES 

Supersonic turbulence in compressible media establishes a complex network of interacting shocks. Converging 
shock fronts locally generate large density enhancements, diverging flows create rarefied voids of low gas density. The 
fluctuations in turbulent velocity fields are highly transient, as the random flow that creates local density enhancements 
can disperse them again. The life time of individual shock-generated clumps corresponds to the time interval for two 
successive shocks to pass through the same location in space, which in turn depends on the length scale of turbulence 
and on the Mach number of the flow. 

The velocity field of turbulence that is driven at large wavelengths is found to be dominated by large-scale shocks 
which are very efficient in sweeping up material, thus creating massive coherent density structures. The shock passing 
time is rather long, and shock-generated clumps can travel quite some distance before begin disrupted. On the 
contrary, when energy is inserted mainly on small scales, the network of interacting shocks is very tightly knit. 
Clumps have low masses and the time interval between two shock fronts passing through the same location is small, 
hence, swept-up gas cannot travel far before being dispersed again. 

The density and velocity structure of three models with large-, intermediate-, and small- wavelength turbulence is 
visualized in Figure ^ It shows cuts through the centers of the simulated volume. As turbulence is stationary, all 
times are equivalent, and the snapshot in the upper panel is taken at some arbitrary time. The lower panel depicts 
the system some time interval later corresponding to 1 /4 shock crossing time through the cube. One clearly notices 
markable differences in the density and velocity field between the three models. 

A. Transport properties in an absolute reference frame 

In order to drive supersonic turbulence and to maintain a given rms Mach number in the flow, we use a random 
Gaussian velocity field with zero mean to 'agitate' the fluid elements at each timestep. However, despite the fact that 
the driving scheme has zero mean, the system is likely to experience a net acceleration and develop an appreciable 
drift velocity, because of the compressibility of the medium. This evolutionary trend is well illustrated in Figure |3 
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FIG. 1: Density and velocity structure of models 2£, 2i, and 2s (from left to right). The panels show cuts through the center of 
the computational volume normal to the three principal axes of the system, after one shock crossing time t^ oss — L/a v « 6.5 
and 1/4 t cross later. Density is scaled logarithmically as indicated in the greyscale key at the upper left side. The maximum 
density is ~ 100, while the mean density is one in the normalized units used. Vectors indicate the velocity field in the plane. 
The rms Mach number is M « 3.1. Large-scale turbulence (2£) is dominated by large coherent density and velocity gradients 
leading a large degree of anisotropy, whereas small-scale turbulence (2s) exhibits noticeable structure only on small scales with 
the overall density structure being relatively homogeneous and isotropic. 



which plots the time evolution of the three components of the mean velocity for models 21, 2i, and 2s, with rms Mach 
numbers M. « 3.1, where turbulence is driven on (a) large (i.e. with small wavenumbers 1 < k < 2), (b) intermediate 
(3 < k < 4), and (c) small scales (with 7 < k < 8). The net acceleration is most pronounced when turbulent energy 
is inserted on the global scales, as in this case larger and more coherent velocity gradients can build up across the 
volume compared to small-scale turbulence. 

The tendency for the zero-mean Gaussian driving mechanisms to induce significant center-of-mass drift velocities 
in highly compressible media can be understood as follows. Suppose the gas is perturbed by one single mode in form 
of a sine wave. If the medium is homogeneous and incompressible, equal amounts of mass would be accelerated in 
the forward as well as in the backward direction. But, if the medium is inhomogeneous, there would be an imbalance 
between the two directions and the result would be a net acceleration of the system. If the density distribution remains 
fixed, this acceleration would be compensated by an equal amount of deceleration after half a period, and the center 
of mass would simply oscillate. However, if the system is highly compressible and the driving field is a superposition 
of plane waves, the density distribution would change continuously (and randomly). Any net acceleration at one 
instance in time would not be completely compensated after some finite time interval later. This will only occur for 
t — > oo assuming crgodicity of the flow. Subsequently, the system is explected to develop a net flow velocity in some 
random direction for t < oo. This effect is most clearly noticeable for long- wavelength turbulence, where density and 
velocity structure is dominated by the coherent large-scale structure. But the effect is small for turbulence that is 
excited on small scales, because in this limit, there is a large number of accelerated 'cells' which in turn compensate 
for another's acceleration. 

The property that the compressible turbulent flows are likely to pick up average drift velocities, even when driven 
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FIG. 2: Time evolution of the mean flow velocity {vi{t))i in models 21, 2i, and 2s. Time t and velocity v are given in normalized 
units (Section IHll . 
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FIG. 3: Time evolution of the diffusion coefficient D(t) = d([fi(t) — fi(0)] 2 )i/dt for models Si, and 2s computed in an 
absolute reference frame. All units are normalized as described in Section UTTl 



by Gaussian fields with zero mean, has implications for the transport coefficients. Figure shows the time evolution 
of the absolute (Eulerian) diffusion coefficients D x , D y , and D z in each spatial direction computed from Equation © 
Note that for stationary turbulence, only time differences are relevant and one is free to chose the initial time. In 
order to improve the statistical significance of the analysis, we obtain D(t) and £f(t) by further averaging over all 
time intervals t that 'fit into' the full timespan of the simulation. 

Due to the (continuous) net acceleration experienced by the system, the quantity &(t) grows faster than linearly 
with time, even for intervals much larger than the correlation time r, i.e. for r <C t < oo. The system resides 
in a superdiffusive regime, where D(t) does not saturate. Instead, D(t) grows continuously with time, which is 
most evident in model 2£ of large-scale turbulence. The ever increasing drift velocity («j(t)), causes strong velocity 
correlations leading continuous growth of the velocity autocorrelation tensor J* %xC(t')dt' . This net motion, however, 
can be corrected for, allowing us to study the dispersion of particles in a reference frame that moves along with the 
average flow velocity of the system. 



B. Transport properties in flow coordinates 

In order to gain further insight into the transport properties of compressible supersonic turbulent flows, we study 
the time evolution of the relative (Lagrangian) diffusion coefficient. In this prescription, 

D'(t) = ^ (6) 

is obtained relative to a frame of reference which comoves with the mean motion of the system («j(t))j on the trajectory 
(n(t))i = Jo(v i (t')) i dt'. Then, g(i - <') = ([(n(t) - (nit)),) - (n(t) - (n(t')) t )} 2 ) t (see Equation©. 

In Figure EJ we show the evolution of D'(t) for each coordinate direction for the complete suite of models. The 
rms Mach numbers range from about 0.5 to 10, each considered for three cases where turbulence is driven on large, 
intermediate, and small scales, respectively. The plots are rescaled such that the time-averaged one-dimensional rms 
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FIG. 4: Time evolution of the diffusion coefficient D'(t) = /dt computed in a reference frame that follows the average flow 
velocity {thick line, axis scaling on the left ordinate), i.e. is centered on (fi(t))i — J*{vi(t')}idt'. Velocity dispersions along 
the three major axes x, y, and z are each normalized to unity using the time-averaged one-dimensional Mach number M (as 
indicated in each plot) together with the given value of the sound speed c s {thin lines, axis scaling on the right ordinate). 
Times are rescaled to the rms shock crossing time through the simulated cube t cross = L/o v — L/(Mc s ). Details for each 
model are given in Table[I] The horizontal gray shaded area indicates the mixing length prediction for t — > oo, and the vertical 
gray and light gray shaded areas show a time interval of r = L/(kA4c s ) and 2r, respectively. For t <g r diffusion should be 
anomalous and coherent, with D'(t) growing linearly with time. The expected behavior from mixing length theory in this 
regime is indicated by the straight line originating at t = 0. It indeed gives a good fit. Note that all models driven on large 
scales (1 < k < 2) exhibit a considerable degree of anisotropy, manifested by different rms Mach numbers M along the three 
principal axes and different values D'(t). For the models with intermediate-scale and small-scale driving anisotropy effects are 
increasingly less important. 
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velocity dispersion a v is normalized to unity (for each direction separately). We also rescale the time t with respect 
to the average shock crossing time scale through the computational volume, t cmss = L/a v . Recall that L = 1, and 
note that a v usually differs between the three spatial directions because of the variance effects, especially in models 
of large-scale turbulence. 

In Figure^] we demonstrate that the magnitude of D'(t) saturates for large time intervals t > t in all directions. In 
a reference frame that follows the mean motion of the flow, diffusion in compressible supersonically turbulent media 
indeed behaves in a normal manner. For small time intervals t < r, however, the system still exhibits an anomalous 
diffusion even with the mean-motion correction. In this regime D'(t) grows roughly linearly with time. For t > r the 
diffusion coefficient D'(t) reaches an asymptotic limit. This result holds for the entire range of Mach numbers studied 
and for turbulence that is maintained by energy input on very different spatial scales. 

From Figure^ we find that diffusion in compressible supersonic turbulent flows follows a universal law. It can be 
obtained by using the rms Mach number (together with the sound speed c s ) as characterizing parameter for rescaling 
the velocity dispersion a v , and the rms shock crossing time scale through the volume i cross = L/(A4c s ) for rescaling 
the time. The normalized diffusion coefficient D'(i) exhibits a universal slope of two at times t < r (i.e. in the 
superdiffusive regime), and approaches a constant value that depends only on the length scale but not on the strength 
(i.e. the resulting Mach number) of the mechanism that drives the turbulence. Even for highly compressible supersonic 
turbulent flows it is possible to find simple scaling relations to characterize the transport properties — analogous to 
the mixing length description of diffusive processes in incompressible subsonic turbulent flows. 

V. A MIXING LENGTH DESCRIPTION 

Incompressible turbulence is often described in terms of a hierarchy of turbulent eddies, where each eddy contains 
multiple eddies of smaller size on the lower levels of the hierarchy, while itself being part of turbulent eddy at larger 
scales [3 Hzl HH • At each level of the hierarchy, an eddy is characterized by a typical lengthscale i and a typical 
velocity v. The typical lifetime of an eddy is its 'turn-over' time r = £/v. This mixing length prescription is an 
attempt to characterize the flow properties in terms of the typical scales £ and v. For example, this classical picture 
defines an effective 'eddy' viscosity /x = p£v, where p is the density. The mixing length £ is interpreted to be the 
turbulent analogue of the mean free path of molecules in the kinetic theory of gases, with v being the characteristic 
velocity of the turbulent fluctuation. 

In such a model, the velocities of gas molecules within an eddy are strongly correlated within a time interval t < r. 
They all follow the eddy rotation; the diffusion process is coherent. However, for t ^> r the velocities of gas molecules 
become uncorrelated, as the eddy has long been destroyed and dispersed. Hence, the velocity autocorrelation function 
vanishes for large time intervals, C(t) — ► for t — > 00. Diffusion becomes incoherent as in Brownian motion or the 
random walk. The diffusion coefficient in the mixing length approach simply is D(t) « 2v 2 t in the regime t < r, which 
follows from replacing f(t) by vt and v(t) by v in Equation [21 As the largest correlation length is the eddy size, f(t) 
is substituted by £ = vt for times t^>r, and the classical mixing length theory yields D(t) s=s 2£v — 2v 2 t = constant. 

Compressible, supersonic, turbulent flows rapidly build up a network of interacting shocks with highly transient 
density and velocity structure. Density fluctuations are generated by locally converging flows, and their lifetimes are 
determined by the time r between two successive shock passages. This time interval is determined by the typical 
shock velocity, which is roughly the rms velocity of the flow, i.e. the Mach number times the sound speed, a v = M.c s . 
It also depends on the length scale at which energy is inserted into the system to maintain the turbulence, which in 
our case is L/k with k being the driving wavenumber and L being the size of the considered region (recall that in our 
models L is unity). This length scale is also the typical traveling distance before two shocks interact with each other. 
As basic ingredients for a supersonic compressible mixing length description we can thus identify: 

shock travel length : £^L/k, (7) 
rms velocity : v ~ a v — Aic s . (8) 

The Lagrangian velocity correlation time scale, r, is analogeous to the time interval during which shock-generated 
density fluctuation remains unperturbed and moves coherently before it is being dispersed by the interaction with a 
new shock front. This time interval is equivalent to the time scale a shock travels along its 'mean free path' £ with an 
rms velocity v. This crossing time is r = £v w a v L/k. For t < r gas molecules can travel coherently within individual 
shock generated density fluctuations, and the diffusion coefficient in the mixing length prescription follows as 

D'{t) w 2v 2 t « 2(j 2 v t . (9) 

D'(t) grows linearly with time with slope 2a 2 . For large times, f>T, D'(t) approaches a constant value, 



D'{t) w 2v 2 t w 2a v L/k 



(10) 
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This mixing length approach (Equations l9land HU|) suggests a unique scaling dependence of the diffusion coefficients 
in supersonic compressible flows on the Mach number M. and on the length scale £ of the most energy containing 
modes with respect to the total size L of the system considered. 

We can use M. (together with the given value of the sound speed) to normalize the rms velocity: a v = Mc s i— > a' v = 1. 
And we can also rescale the time with respect to the rms shock crossing time scale through the total volume, which 
is t CIOSS — L/a v — L/(Aic s ) = t soun d/A4 with i SO und = L/c s being the sound crossing time, so that 1 1— > t' = t/t CIOSS . 
From this normalization procedure, we get D'(t) i— > D"(t') — D'(t) M.c s L and obtain the following universal profile 
for the diffusion coefficient, 



D"(t') = 2t' for 
D"(t') = 2/k for 



t' < 1/fc and 
t' > 1/fc . 



(11) 
(12) 



Note that this result holds for each velocity component separately, as the results in Figure 0] indicate. In this 
case o~ v stands for a x , a y , or a z in Equations 151 and I1UI and it hold for the total diffusion coefficient, when using 



l 1 / 2 instead. 
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FIG. 5: Comparison between mixing length predictions and numerical models. At the left we plot he normalized, mean-motion 
corrected diffusion coefficient D"(t') for t — > oo, and at the right its slope dD" {t') / dt' for t' <C 1/fc. For each suite of models, 
large-scale, intermediate-scale, and small-scale turbulence, respectively (as indicated by the forcing wavenumber fc at the top 
of each plot), we separately show the three velocity components (as indicated at the bottom). The different Mach numbers in 
each model suite are denoted by different symbols (as identified at the right-hand side of each plot). The dotted lines give the 
corresponding prediction of the mixing- length theory, D"(t') — 2/k and dD" (t 1 ) /dt' — 2, respectively, where we take fc to be 
the maximum wavenumber of the forcing scheme (indicated at the top of each plot). 



The validity of the mixing length approximation is quantified in Figure [S] which plots the mixing length predictions 
against the values obtained from the numerical models. For large-scale and intermediate scale turbulence, the mixing 
length approach gives very satisfying results, only for small-scale turbulence it underestimates the diffusion strength. 
This disparity probably has to do with the numerical resolution of the code, in the sense that driving wavenumbers of 
fc w 8 come close to the dissipation scale of the method and hence the inertial range of turbulence is limited 23] . That 
limitation leads to an effective driving for the turbulent motion on somewhat larger scales than 1/8. Consequently, it 
leads to a stronger diffusion, i.e. somewhat larger diffusion coefficients than those predicted by Equations 1111 and 1121 
The same numerical effects also account for the slightly shallower slope of D'(t) for t <C r for models 7 < fc < 8. 

Figures 01 and |S] indicate that the classical mixing length theory can be extended from incompressible (subsonic) 
turbulence into the regime of supersonic turbulence of highly compressible media. In this case, driving length £ and 
rms velocity dispersion a v = A4c s act as characteristic length and velocity scales in the mixing length approach. Note, 
that this only applies to mean-motion corrected transport. In general (i.e. in an absolute reference frame), supersonic 
turbulence in compressible media leads to superdiffusion as visualized in Figure 
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VI. SUMMARY 

We studied diffusion processes in supersonically turbulent, compressible media. To drive turbulence and maintain 
the desired rms Mach number in the flow, we insert energy into the system at a pre-specified rate and over a given 
spatial scale using random Gaussian fields. In our numerical models, the adopted magnitude of the rms Mach numbers 
range from A4 = 0.5 to M. = 10, and turbulence was driven on large, intermediate, and small scales, respectively. 

Supersonic turbulence in compressible media establishes a complex network of interacting shocks. Converging shock 
fronts locally generate large density enhancements, diverging flows create voids of low gas density. The fluctuations in 
turbulent velocity fields are highly transient, as the random flow that creates local density enhancements can disperse 
them again. 

Due to compressibility, supersonically turbulent flows will usually develop noticeable drift velocities, especially when 
turbulence is driven on large scales, even when it is excited with Gaussian fields with zero mean. This tendency has 
consequences for the transport properties in an absolute reference frame. The flow exhibits super-diffusive behavior 
(see also Q). However, when the diffusion process is analyzed in a comoving coordinate system, i.e. when the induced 
bulk motion is being corrected, the system exhibits normal behavior. The diffusion coefficient D(t) saturates for large 
time intervals, t — > oo. 

By extending classical mixing length theory into the supersonic regime we propose a simple description for the 
diffusion coefficient based on the rms velocity v of the flow and the typical shock travel distance £, 

D'(t) = 2v 2 t for t<s:£/v, 
D'(t) = 2v£ for f>l/v. 

This functional form may be used in those numerical models where knowledge of the mixing properties of turbulent 
supersonic flows is required, but where these flows cannot be adequately resolved. This is the case, for example, in 
astrophysical simulations of galaxy formation and evolution, where the chemical enrichment of the interstellar gas and 
the distribution and spreading of heavy elements produced from massive stars throughout galactic disks needs to be 
treated without being able to follow the turbulent motion of interstellar gas on small enough scales relevant to star 
formation 0, l50| . Our results furthermore are directly relevant for understanding the properties of individual star- 
forming interstellar gas clouds within the disk of our Milky Way. These are dominated by supersonic turbulent motions 
which can provide support against gravitational collapse on global scales, while at the same time produce localized 
density enhancements that allow for collapse, and thus stellar birth, on small scales. The efficiency and timescale of 
star formation in galactic gas clouds depend on the intricate interplay between their internal gravitational attraction 
and their turbulent energy content |30j. The same is true for the statistical properties of the resulting star clusters. 
For example, the element abundances in young stellar clusters are found to be very homogeneous [59|. implying that 
the gas out of which these stars formed must be have been chemically well mixed initially. On the basis of the results 
discussed here, this observation can be used to constrain astrophysical models of interstellar turbulence in star-forming 
regions. Understanding transport processes and element mixing in supersonic turbulent flows thus is a prerequisite 
for gaining deeper insight into the star formation phenomenon in our Galaxy. 
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